A method of blasted rock image segmentation based on improved watershed algorithm

It is of great theoretical significance and practical value to establish a fast and accurate detection method for particle size of rock fragmentation. This study introduces the Phansalkar binarization method, proposes the watershed seed point marking method based on the solidity of rock block contour, and forms an adaptive watershed segmentation algorithm for blasted rock piles images based on rock block shape, which is to better solve the problem of incorrect segmentation caused by adhesion, stacking and blurred edges in blasted rock images. The algorithm first obtains the binary image after image pre-processing and performs distance transformation; then by selecting the appropriate gray threshold, the adherent part of the distance transformation image, i.e., the adherent rock blocks in the blasted rock image, is segmented and the seed points are marked based on the solidity of the contour calculated by contour detection; finally, the watershed algorithm is used to segment. The area cumulative distribution curve of the segmentation result is highly consistent with the manual segmentation, and the segmentation accuracy was above 95.65% for both limestone and granite for rock blocks with area over 100 cm2, indicating that the algorithm can accurately perform seed point marking and watershed segmentation for blasted rock image, and effectively reduce the possibility of incorrect segmentation. The method provides a new idea for particle segmentation in other fields, which has good application and promotion value.

Blasting is widely used in mining and civil engineering due to its economy and efficiency [1][2][3][4][5] . As an important technical indicator of blasting effectiveness, blasted block size distribution directly affects the cost and efficiency of subsequent shoveling, crushing and grinding processes, and also provides a necessary basis for blasting parameter optimization [6][7][8][9][10] . Therefore, it is of theoretical significance and practical value to establish a fast and accurate detection method for particle size of rock fragmentation to guide blasting construction and improve blasting efficiency. Blasted rock piles are characterized by large scale, serious adhesion and irregularly shaped rock clumps, large differences in particle size, and small differences in grayness, which make it difficult to accurately measure the particle size of blasted rocks 11,12 . The existing measurement methods can be summarized into two categories: three-dimensional (3D) point cloud data segmentation measurement and two-dimensional (2D) image segmentation measurement. The 3D point cloud data is mainly obtained using 3D laser scans or a large number of high-resolution digital photos taken by Unmanned Aerial Vehicle (UAV) 13,14 . Han and Song 15 used stereophotogrammetry for 3D modeling of blasted rock piles to obtain surface block dimensions and corrected for errors in rock fast due to stacking by indoor tests, and validated the applicability of the method in small-scale tests in the field. Although the creation of point cloud data using stereophotogrammetry is more effective and cheaper, the large number of images received takes a lot of time to process and convert into point cloud data 16 . 3D laser scanners are widely used in mine surveying due to its directly, quickly and captures 3D geometry in detail [17][18][19][20][21] . Engin et al. 22 used a 3D laser scanner to obtain a 3D view of a rock piles of about 13 cm and used morphological methods to determine the position of the rock block and to correct the surface of the rock block, and finally used nonlinear order statistical filtering and histogram analysis to determine the blasted block size distribution of the rock piles, and by comparing the results of image analysis, the results obtained using this method were proved to be sufficiently reliable and accurate. Wang et al. 23 used 3D laser scanning technology to obtain blasted rock piles point cloud data, and used the Voxel Cloud Connectivity Segmentation algorithm improved by discrete features to solve the influence of point clouds on the surface of small particles of ore on block recognition, and used the Locally Convex Connected Patches algorithm improved by plane fitting to solve the problem of over-segmentation of large rock blocks, and verified the generality and accuracy of the method

Image acquisition
Limestone quarries and granite quarries located in Huizhou, China were selected for the study, are both operated by Guangdong Xiyuan Blasting Technology Co., Ltd. Therefore, the blasting design parameters are the same, except for the charge quantity. The step height of the quarry is about 12, the diameter of the gun hole is 140 mm, using a digital electronic detonator to hole-by-hole detonation. Figures 1 and 2 show the 3D images established by aerial photography of a limestone quarry and a granite quarry after a certain blast using a UAV, respectively.
The UAV used the DJI Phantom 4 RTK, which is equipped with a multi-frequency and multi-system RTK module. The camera is synchronized with the RTK module µs-level time to provide real-time cm-level positioning data without the need to deploy image control points. The UAV is equipped with a high-precision anti-shake tripod head camera that supports up to 20 megapixels of still photo shooting. To avoid errors caused by rock block overlap, image acquisition of blasted rock piles is performed using tilt photogrammetry, where the camera orientation is perpendicular to the blasted rock piles surface. Figure 3 shows the flight schematic of the UAV. Firstly, determined point A, B and C. Then, the outward expansion was carried out from point C to A and B respectively until the flight range completely covered the rock piles, and the UAV automatically planned the flight route according to the overlap rate. To ensure clarity in images, the flight altitude was set to 25, which is the minimum flight altitude of the UAV. The collected images of limestone and granite blasting blocks are shown in Fig. 4

Proposed methods
Image pre-processing. Due to the more complex and dustier environment of the quarry, there is more serious noise in the images of blasted rock piles collected through the camera, and the rock piles are heavily stacked and adhered to each other, with small background difference degree and inconspicuous color information. To effectively segment the rock blocks, pre-processing is required for the blast rock piles images. There are many methods of image pre-processing 47,48 , and this study adopts the more commonly used methods in the field of rock segmentation. Firstly, the grayscale image is de-noising by bilateral filtering. Secondly, the rock block     Figure 5 shows the effect of limestone blasting rock block image after pre-processing. The Phansalkar method 49 is an adaptive local thresholding method based on local image properties, which processes each pixel (x, y) in an image by considering a w × w window with the pixel as the center on that window with a grayscale mean is m(x, y) and standard deviation is s(x, y) , then the local threshold T(x, y) for the pixel is: where p and q are constants. The Phansalkar method is flexible in that it determines the selection of the threshold value based on the magnitude of the local mean and the standard deviation, and by adjusting the values of the parameters p , q and k , different processing results are obtained. the processing effects are shown in Fig. 5c.
To show indicate the effectiveness, the method is compared with the Otsu method 50 , and the results of the processing of Fig. 5b using the Otsu method are given in Fig. 6, from which it can be seen that the Otsu method produces three kinds of errors: (1) incorrectly dividing regions into backgrounds, such as darker regions and zones obscured by shadows; (2) incorrectly dividing small rocks into backgrounds; (3) dividing part of the inner regions of large rocks into backgrounds. Some of the error areas are shown in the red box part of the image, which can seriously affect the subsequent rock identification. The Phansalkar method can accurately distinguish the background from the rock particles with satisfactory results. Figure 5d shows the binary image of the blasted rock piles after morphological optimization. Comparing this image with Fig. 5c, we can see that morphological optimization can eliminate the small "holes" and noise points in the binary image and smooth the target edges, but some of them still cannot be completely removed. Figure 5e uses the area filtering elimination method to achieve more satisfactory results, but there are still some noise points, which affect the rock segmentation. When the area threshold increases, it causes some edges to be removed. As in the red box in Fig. 5d, small rock block contour are eliminated, which has an impact on the segmentation results.
Principle of watershed segmentation algorithm based on distance transformation. The principle of the watershed algorithm is to visualize an image as a 3D topographic image 34 . In the terminology of "topography", three types of points are considered. In Fig. 7: ① a local min point (min value surface), which corresponds to the lowest point of a basin; ② points at other locations in the basin; ③ points at the edge of the basin, where the basin meets other basins. For a specific regional minimum, the set of points meeting condition ② is called the catchment basin or watershed for that minimum, and the points meeting condition ③ form the front line of the ground called the division line or watershed line. The main goal of the algorithm is to find the watershed line, which is the contour of the rock block image. Assuming that a hole is punched at the minimum of each area, we let the water pass through the hole to flood the entire terrain at a uniform rate. As the water rising in the different catchment basins aggregates, a dam is built to stop this aggregation until the water floods the highest point of the topographic image, and the boundaries of these dams are the dividing lines of the watershed.
The traditional watershed algorithm is usually marker less segmentation, and the input object is a gradient image, which is based on the luminance change, and it only reflects the edge information of the image, which will result in unreasonable segmentation because of its noise-sensitive feature, leading to the following disadvantages of the watershed algorithm: ① the noise in the original image causes the segmentation contour shift; ② the image with low contrast, the contour of the target region is easily lost when segmenting; ③ there are many meaningless www.nature.com/scientificreports/ local minima in the image. Therefore, the most common practice is to perform a distance transformation on the binary image. The distance from each pixel in a binary image to its nearest zero-valued pixel is called the distance transform 51 . Suppose a binary image I with a target set O and a background set B and a distance image D. The distance transformation is defined as in Eq. (2).
The Euclidean distance is generally chosen as the disf () . The calculation method is as follows: www.nature.com/scientificreports/ Defects of watershed segmentation algorithm based on distance transformation. Figure 8 gives the distance transformed image of the binary image of the blasted rock piles. From Fig. 8b, the distance image of the binary image is similar to image skeletonization and still retains the general shape. Figure 9 shows the distance transformation detail image and seed point image of a rock block, where (a) is a rock block in the image; (b) is the distance transformation image; (c) is the local extreme value point (seed point). From Fig. 9c, it can be seen that there are multiple extreme value points, i.e., redundant seed points, in   www.nature.com/scientificreports/ the center of the block (in the red box), which will cause more serious over-segmentation, and there are also more extreme value points in the sticky part of the two blocks (in the yellow box), which seriously affects the segmentation effect of the adhesion block. For this problem, the commonly used methods are mainly to merge adjacent maximal points and merge expanded maximal points. But as shown in Fig. 9c, the distance between the maximal points inside the rock block and the maximal points in the adhering part of the rock block is large. The above method is difficult to achieve accurate merging. Moreover, for the complex image, the distribution of the maximal points is irregular, and the above maximum point elimination method is less effective.

Solidity of rock block contour.
In the segmentation process of some special adherent objects, such as adherent cells, spherical particles, etc., some prior knowledge, such as the color of cell nuclei, shape of spherical particles and others, can be used in turn to correct the seed points. However, for the segmentation of blasted rock block, there is no a priori knowledge available other than the shape, and the shape of blasted rock block is irregular polygons. To investigate the rule of blasted rock block shape, this study performs manual segmentation of blasted rock images, as shown in Fig. 10.
Statistical information of rock blocks, including area, contour convex hull, solidity of rock block, through image contour detection technology, and it should be noted that the solidity is the ratio of area to the contour convex hull. Due to the limitation of space and considering the large segmentation error of small rocks, only part of the information of larger blocks is shown in Tables     www.nature.com/scientificreports/ as a seed point if the contour solidity is greater than the set solidity threshold. As shown in Fig. 13, (a) is a 3D schematic image of the distance image, and the height is the gray value of the corresponding coordinate; (b1) is a 3D schematic image of the distance image when the gray threshold value is 50, and (b2) is its corresponding 2D cross-section; (c1) is a 3D schematic image of the distance image when the gray threshold value is 100, and (c2) is its corresponding 2D cross-section. When the grayness threshold is too small, although smaller blocks can be split and seed point marked, more large rock blocks are not divided, as can be seen in Fig. 13b1,b2. And most of the contour solidity is less than the solidity threshold value, which does not work as seed points. It can be seen from Fig. 13c1,c2 that when the gray threshold is too large, only the large size rock blocks can be segmented. Due to the large size difference, a single gray threshold does not satisfy the seed point marking of rock blocks, so multiple gray thresholds need to be processed and the solidity is calculated for the contours in the 2D image after each gray threshold is processed. If it is greater than the solidity threshold, the interior of the contour is filled and marked as a seed point.
It should be noted that since this seed point marking method is based on the distance transform image, there may be two cases as follows: ① When the solidity of the contour of a background region in a binary image is too high, a depressed structure as in Fig. 14a will be formed in this background region after the distance transformation process, when the contour is taken using the gray threshold, a contour as in Fig. 14b will be formed, and the solidity of this contour is greater than the solidity threshold, which will cause wrong segmentation if the contour is marked as a seed point. ② As shown in the red marker in Fig. 15, this noise is caused by the wrong differentiation of Phansalkar. Obviously, morphological optimization and area filtering do not eliminate this noise, and the contour is smaller than the solidity threshold value. Therefore, when determining whether the contour is a seed point or not, it is necessary to determine whether there is a background area or noise within the contour, and if it exists, it is necessary to determine whether the contour is generated by the background area or noise.
The pseudo code of the seed point marking method based on the solidity of rock block contour is shown in steps 8-12 of algorithm 1.
The adaptive watershed algorithm based on the solidity of rock block contours. Achieving accurate segmentation of blasted rock images requires a complete process, which mainly consists of image denoising, histogram equalization, image binarization, morphological optimization, distance transform, seed point marking and watershed segmentation. Based on the above research, an adaptive watershed segmentation algorithm is proposed for blasted rock piles images based on the solidity of rock contours. The method performs adaptive segmentation based on the gray feature of blasted rock piles image and the rock contour features      Fig. 4 as the test objects, the applicability and accuracy of the adaptive watershed algorithm based on rock block contour solidity for blasted rock images are verified. Firstly, the images are preprocessed, and the binary images of limestone blasted rock is shown in Fig. 8a, and the binary images of granite blasted rock are shown in Fig. 16. Then the blasted rock images are segmented using the method proposed in this study. Fig-Figure 16. Granite blasted rock binary image.  Figs. 17, 18 and 19, it can be seen that the seed point marking method based on the solidity of rock block contours can effectively mark the rock blocks, especially the larger ones. Comparing the seed point image of the blasted rock piles with its pre-processed binary image, it is found that the seed point marking method can avoid the effect of noise inside the binary image on the segmentation, as shown in the red box in Fig. 16. However, the method also has some problems.

Analysis of seed point marking results. From the seed point images in
(1) Some of the severely adhered blocks are not separated, as shown in the red boxes in the seed point images of Figs. 17, 18 and 19. It can be seen from the image that the problem is likely to arise when a large block heavily adheres to a small block and the difference in grayscale is small, due to the result that there are no local minima in the small block after the distance transformation and that solidity of the contours formed by the small and large rock blocks during seed point marking is greater than the solidity threshold.  www.nature.com/scientificreports/ By comparing the seed point images with different solidity thresholds, it can be found that with the increase of solidity threshold, the case of no segmentation of the adherent rock block gradually decreases, but the case of multiple seed point marking of the rock block gradually increases, so it is necessary to select a suitable solidity threshold for seed point marking of the rock block to achieve the best segmentation effect.

Analysis of image segmentation results of blasted rock piles. From the segmentation results in
Figs. 17, 18 and 19, it can be seen that the adaptive watershed algorithm based on the solidity of rock block contours can achieve more accurate segmentation of the severely stacked and adhered blast rock images, especially the limestone blast rock images with less noise. In contrast, the segmentation effect of granite blast rock images is poor compared with that of limestone blast rock images due to the influence of rock block shadows and other problems. In order to evaluate the segmentation results with quantitative indexes, this study uses the manual segmentation image as the segmentation criterion to evaluate the segmentation effect of two blasted rock images, in which the manual segmentation image of limestone blasted rock is shown in Fig. 10a and the manual segmentation image of granite blasted rock is shown in Fig. 10b. Figure 20 show the cumulative distribution curves of the area of rock blocks in the images of limestone and granite blasted rock piles. The calculation formula is shown in Eq. (4).
where P is the cumulative area ratio of rock blocks; S total is the sum of the areas of all identified rock blocks across the blasted rock image; and S x is the cumulative area of rock blocks in part x of area classification. www.nature.com/scientificreports/ As can be seen from Fig. 20, for rock segmentation of limestone blasted rock piles images, the agreement between the cumulative distribution curves of the area of the algorithm proposed in this study and the manual segmentation results gradually increases with the increase of the solidity threshold. but rock segmentation of granite blasted rock images becomes progressively worse as the solidity threshold increases. From Fig. 20c1, it can be seen that there are some differences in the segmentation results of limestone rocks with area ranges between 3600-4300 and 5500-6000 cm 2 . From Fig. 20a2, it can be seen that there are some differences in the segmentation results of granite rocks in the area range between 4000-5000 and 6000-7000 cm 2 , and the maximum area rock segmentation error is larger. Table 3 gives a comparison of the three characteristic area parameters, Area 20 , Area 50 and Area 80 , at different solidity thresholds, which represent the areas corresponding to cumulative area distribution percentages of 20%, 50%, and 80% respectively. As can be seen from the table, like the solidity threshold increases, the Area 20 error for limestone gradually increases, while the Area 50 and Area 80 errors gradually decrease. The opposite is true for granite, where the Area 20 error decreases with increasing solidity threshold, while the Area 50 and Area 80 errors gradually increase. Tables 4 and 5 give the comparison of the number of rocks in different area zones for the image segmentation of limestone and granite blasted rock piles, respectively, and it should be noted that only the number of rocks with an area of 100 cm 2 or more is counted. As can be seen from Table 4, the accuracy of block segmentation for limestone blasted rock images above 100 cm 2 is above 95.80%. with the increase of solidity threshold, the number of rocks in the 100-1000 cm 2 interval gradually increases, and the rocks above 1000 cm 2 are less affected by the solidity threshold. Table 5 shows that the accuracy of rock segmentation for granite blasted block images above 100 cm 2 is above 95.65%, and the number of rocks in the 100-1000 cm 2 interval gradually increases with the increase of solidity threshold, while the number of rocks above 1000 cm 2 gradually decreases. This is mainly due to the fact that as the solidity threshold increases, the area of large rock seed points gradually decreases and even splits into multiple seed points, resulting in an increase in the number of small rocks and a decrease or no  Therefore, among the three solidity thresholds, the best solidity threshold interval for limestone blasted rock image segmentation is 0.85-0.9, the best solidity threshold interval for granite blasted rock image segmentation is 0.8-0.85, considering the cumulative area distribution curve, characteristic area parameters and the number of rock blocks. In general, the algorithm proposed in this paper is consistent with the manual segmentation results under the condition that a suitable solidity threshold is selected.
Comparison with current methods. To evaluate the performance of the algorithm proposed in this study, the algorithm was compared with two other watershed improvement methods for rock segmentation, as described in the literature 45,46 , respectively. Figures 21 and 22 show the segmentation results using the segmentation methods from the literature 45,46 , respectively. Figures 23 and 24 show the cumulative area distribution curve of the literature 45,46 segmentation results, respectively. The three characteristic area parameters, Area 20 , Area 50 and Area 80 , are shown in Table 3. The comparison of the number of rocks in different area for the segmentation of limestone and granite blasted rock images is shown in Tables 4 and 5, respectively.
As can be seen from Figs. 23a and 24a, the method of literature 46 gives better results for limestone compared to the method of literature 45 , with cumulative area distribution curve over 1000 cm 2 is more similar to the result of manual segmentation. As can be seen from Figs. 23b and 24b, the results of the method of literature 45 are better than those of literature 46 for granite. From Tables 3, 4 and 5, it can be seen that among the three segmentation methods, the method proposed in this study has the best segmentation results, and its segmentation results will be closer to the manual segmentation after selecting a suitable solidity threshold.

Conclusion
In this study, to avoid errors caused by rock block overlap, image acquisition of blasted rock piles is performed using tilt photogrammetry. Then, the accurate rock target area was obtained through the image pre-processing process. Finally, the binary image is segmented using the adaptive watershed segmentation algorithm of blasted rock image based on rock block shape. This method enables automatic segmentation of blasted rock piles rock     www.nature.com/scientificreports/ blocks is achieved. We compared the performance of the proposed method with the manual sieving results and with three methods in the literature 45,46 . The main conclusions are as follows: (1) The pre-processing process of blasted rock image and its influence on the segmentation results are described, and the Phansalkar method is introduced. By comparing the results with the Otsu method, it is proved that the method is more applicable to the binarization of blasted rock images. (2) The principle of the watershed algorithm based on distance transformation is described, and the reasons why the algorithm is prone to serious over-segmentation are analyzed. The study revealed a high-solidity of the blasted rock block contour by analyzing the shape of blocks, a seed point marking method based on the solidity of rock block contours is proposed, and forms an adaptive watershed segmentation algorithm for blasted rock images based on block shapes. The algorithm solves the problem of mis-segmentation of blasted rock images caused by severe adhesion and large differences in particle size. (3) The algorithm can effectively mark the seed points of blasted rock blocks and avoid the effect caused by noise inside the rock blocks of binary images, the segmentation results are highly similar to the area cumulative distribution curve of the manual segmentation results, and the errors of the different indexes are smaller than the segmentation methods adopt in the literature 45,46 , which proves the effectiveness and superiority of the algorithm for rock block segmentation of blasted rock piles images.
The adaptive watershed algorithm based on the solidity of rock block contours proposed in this study considers the contour properties of rock blocks in the segmentation, resulting in a more significant improvement in segmentation accuracy. However, it should be acknowledged that the method requires different solidity thresholds for different block image segmentation to achieve the best segmentation results. In future research, solidity thresholds for more kinds of rock will be investigated and determined, making the proposed method more widely applied.